#include <stdio.h>
#include "spline_lib.h"
#include "randomlib.h"
#include "step_down.h"
#include <unistd.h>

double fun(double x, int i);
void dat_generate(struct L2_1D_DATA*, int, int, double);
int par_args(int ac, char** av);

/* Global variables for argument passing */
static int inform;
static int nobs;
static int num_fn;
static char method;


int main(int ac, char** av)
{

  int i;
  struct L2_1D_DATA *  dat_actual;
  par_args(ac, av);
  /* Initialize the data_actural. */
  dat_actual = data_1d_initialize( nobs, &inform);
  dat_generate(dat_actual, nobs, num_fn, 0);

  for(i = 0; i < nobs; i++)
    printf("%6f %6f\n", dat_actual->xdata[i], dat_actual->ydata[i]);

  return 0;
}

int par_args(int ac, char** av)
{
  int c;

  while( (c = getopt (ac, av, "n:f:m:")) != -1 )
    {
    switch (c)
      {
      case 'n': 
	nobs = atoi(optarg);  break;
      case 'f': 
	num_fn = atoi(optarg);  break;
      case 'm': 
	method = optarg[0];  break;
      default:
	exit(1);
      }
    }
  return 0;
}

void dat_generate(struct L2_1D_DATA* dat, int num_data, int num_fn, double sigma){

  int i;
  double x;
  double h = 1./(dat->nobs - 1);

  int seed1 = 100 + 10*num_data;
  int seed2 = 1000 + 5*num_data;

  RandomInitialise(seed1, seed2);

  for(i=0; i<dat->nobs; i++){
    dat->weight[i] = 1;
    dat->xdata[i] = i*h;
    x = dat->xdata[i];
    dat->ydata[i] = fun(x, num_fn) + sigma*RandomGaussian(0., 1.);
  }
}

double fun(double x, int i){
	double y= 0.;
	switch(i){
		case 1: y = 10*(sin(x*x*x*12  ));  break;
		case 2: y = 1/(.1 + sin(x*x*8 )*sin(x*x*8 )); break;
		case 3: y = fabs(cos(x*x*6)); break;
		case 4: y = 3*exp( -50*(x-.1)*(x-.1)) + 6*exp( -100*(x-.7)*(x-.7));break;
		case 5: y = 10*(pow(fabs(x-.2), 3) - pow(fabs(x-.3), 3) + 4*pow(fabs(x-.5), 3) - pow(fabs(x-.8), 3)); break;
		case 6: y = (fabs(x-.2) + fabs(x-.3) + fabs(x-.5) + fabs(x-.8)); break;
		case 7: y = sin((1./(x*x+.1)) + 1./(x-1.2)); break;
		case 8: y = 8*exp( -100*(x-.1)*(x-.1)) + 
					6*exp( -200*(x-.6)*(x-.6))+ 
					10*exp( -150*(x-.9)*(x-.9));break;
		case 9: y = 2*sin(8*exp( -100*(x-.1)*(x-.1))) + 
					6*exp( -200*(x-.6)*(x-.6))+ 
					10*exp( -150*(x-.9)*(x-.9));break;
	
	}
	return y;
}
